/*********************************************************************************************************
*							
* VARIABLE NAMES
* newid1			-> state, 1=Gujarat, 3=Uttar Pradesh
* newid2			-> district
* village2			-> village2 (numeric), village2char (character)
* detailmodel 		-> detailing status, 0=no, 1=yes
* timepoint 		-> 0=midline, 1=endline
* childgendermodel 	-> 0=female, 1=male
* childagemodel 	-> continuous 0-4
* childnewagemodel  -> categorical
* gendermodel 		-> gender of caregiver, 0=female, 1=male
* agemodel 			-> age of caregiver, 1=18-24, 2=25-29, 3=30-39, 4=40+
* educationhhmodel 	-> education level of caregiver, 0=never been to school, 1=ever been to school
* media 			-> whether caregiver recognizes a specific media message, 0=no, 1=yes
* orsmodel 			-> MODEL 1 OUTCOME, whether the child received any ORS, 0=no, 1=yes
* orszincmodel 		-> MODEL 2 OUTCOME, whether the child received both ORS and zinc, 0=no, 1=yes
* urbanmodel 		-> urban or rural, 0=rural, 1=urban;
* wealthmodel 		-> wealth status, 1=lowest, 2=second, 3=middle, 4=fourth, 5=highest
* 
**********************************************************************************************************/

/* NOTE Make sure the filepath is correct for your particular file location of the dataset */

/* proc contents will display a list of all the variables in this dataset */
proc contents data='h:\chai\midendchildcombined2'; run;

/*********************************************************************************
* table x                                                                        *
* descriptive characteristics stratified by state, time point, and detail status *
*********************************************************************************/

/* proc freq will give you n and % of a single variable. restrict the dataset to subgroups with a 'where' statement */
proc freq data='h:\chai\midendchildcombined2';
table timepoint orsmodel;
*where (timepoint=0) and (newid1=3);
run;

/* two variables joined by * will create a cross-tab */
proc freq data='h:\chai\midendchildcombined2';
table orszincmodel*detailmodel / norow nopercent chisq;
where (timepoint=0) and (newid1=3);
run;

/* surveyfreq accomodates complex survey designs - the cluster statement below accounts for clustering at the village level */
proc surveyfreq data='h:\chai\midendchildcombined2';
table orszincmodel*detailmodel / chisq;
where (timepoint=0) and (newid1=3);
cluster village2char;
run;

/**************************************************************
* table 2                                                     *
* descriptive characteristics stratified by state and outcome *
**************************************************************/

/* variables by ORS/ORSzinc status. use proc freq to fill in columns for row % and n */
proc freq data='h:\chai\midendchildcombined2';
table source*orszincmodel / nopercent nocol;
where newid1=3 ;
run;

/* glimmix bivariate logistic. this produces bivariate odds ratios, 95% CI, and pvalues */
/* for each bivariate test, change the covariate in the class statement (mind the reference group), the covariate and 
   the outcome in the model statement, and the state in the where statement*/
proc glimmix data="h:\chai\midendchildcombined2" method=laplace pconv=0.0001; 
	class	source (ref='1') newid2 village2char;

	model orszincmodel (desc) =	source / dist=bin link=logit solution oddsratio cl;

	random intercept / subject = newid2(village2char);

	where newid1=3;
run;

/************************************************************************************************
* table 3                                                                                       *
* Model 1 = main effects + covariates                                                           *
* Model 2 = main effects + covariates + time*detail interaction                                 *
* Model 3 = main effects + covariates + time*detail interaction + time*detail*media interaction *
************************************************************************************************/

/* there are four possible combinations of outcome/state for each model
   ORS Gujarat, ORS Uttar Pradesh, ORSzinc Gujarat, ORSzinc Uttar Pradesh
   all you need to change is the outcome in the model statement, and the state in the where statement */

*Model 1 main effects + covariates;
proc glimmix data="h:\chai\midendchildcombined2" method=laplace pconv=0.0001; 

	class	timepoint (ref='0') detailmodel (ref='0') media (ref='0') newid2 village2char 
			childgendermodel childnewagemodel educationhhmodel urbanmodel wealthmodel source; 

	model orsmodel (desc) =		timepoint detailmodel media
								childgendermodel childnewagemodel educationhhmodel urbanmodel wealthmodel source
								/ dist=bin link=logit solution oddsratio cl;

	random intercept / subject = newid2(village2char);

	where newid1=3;

run;

*Model 2 main effects + covariates + 2 way interaction;
*NOTE SAS does not generate odds ratios for interaction terms. To express the overall effect of an interaction term in odds,
 exponentiate the coefficient estimate in the 'solution for fixed effects' table by hand;
proc glimmix data="h:\chai\midendchildcombined2" method=laplace pconv=0.0001; 

	class	timepoint (ref='0') detailmodel (ref='0') media (ref='0') newid2 village2char
			childgendermodel childnewagemodel educationhhmodel urbanmodel wealthmodel source; 			

	model orsmodel (event='1') =	timepoint detailmodel media
										timepoint*detailmodel 
										childgendermodel childnewagemodel educationhhmodel urbanmodel wealthmodel source
										/ dist=bin link=logit solution oddsratio cl; 

	*lsmeans timepoint*detailmodel / slicediff=detailmodel oddsratio cl;

	/* NOTE we are not reporting any results from the lsmeans statement, but it provides additional estimates for 
	   variables involved in interaction terms at different levels of other interaction variables */
	
	random intercept / subject = newid2(village2char);

	where newid1=3;

run;



*Model 3 main effects + 2way + 3way;
*NOTE SAS does not generate odds ratios for interaction terms. To express the overall effect of an interaction term in odds,
 exponentiate the coefficient estimate in the 'solution for fixed effects' table by hand;
*NOTE Additionally, because media=0 for all observations at time=0, SAS also does not generate odds ratios for any single variable
 that is involved in an interaction that includes both time and media (i.e., time*detail*media). Similar to the process above, 
 solve for the odds of time, detail, media, and the interactions, by exponentiating the coefficient estimates in the
 'solutions for fixed effects' table by hand;
proc glimmix data="h:\chai\midendchildcombined2" method=laplace pconv=0.0001; 

	class	timepoint (ref='0') detailmodel (ref='0') media (ref='0') newid2 village2char
			childgendermodel childnewagemodel educationhhmodel urbanmodel wealthmodel source; 			

	model orszincmodel (event='1') =	timepoint detailmodel media 
										timepoint*detailmodel
										timepoint*detailmodel*media
										childgendermodel childnewagemodel educationhhmodel urbanmodel wealthmodel source
										/ dist=bin link=logit solution oddsratio cl; 

	*lsmeans timepoint*detailmodel*media / slicediff=detailmodel oddsratio cl; 
	*lsmeans timepoint*detailmodel*media / slicediff=timepoint oddsratio cl; 
	*lsmeans timepoint*detailmodel*media / slicediff=media oddsratio cl; 

	random intercept / subject = newid2(village2char);
	where newid1=1;

run;

/******
* END *
*******




















